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Abstract 

The properties of incompressible, unsteady, interactive, boundary layers are examined for the model hyper- 
sonic boundary layer of Brown et al (1974,1975) and internal flow past humps or, equivalently, external flow 
past short-scaled humps. Using a linear high frequency analysis, it is shown that the domains of dependence 
within the viscous sublayer may be a strong function of position within the sublayer and may be strongly in- 
fluenced by the pressure displacement interaction, or the prescribed displacement condition. Detailed calcula- 
tions are presented for the hypersonic boundary layer. This effect is found to carry over directly to the fully 
viscous problem as well as the nonlinear problem. In the fully viscous problem, the non-local character of 
the domains of dependence manifests itself in the sub-characteristics. Potential implications of the domain 
of dependence structure on finite difference computations of unsteady boundary layers are briefly discussed. 

1. Introduction 

The phenomena of unsteady separation and vortex eruption appear to involve a number of possible high 
Reynolds number structures, including: classical boundary layers and singularities (Van-Dommelen & Shen 
(1980), Smith Si Burggraf (1985)), interactive boundary layers and interactive singularities (Brotherton-Rat- 
diffe Si Smith (1987), Smith (1988)), and the introduction of inviscid effects and normal pressure gradients 
(Elliott et al (1983)). As noted above, at some point the vortex eruption may pass through a stage in which 
viscous-inviscid interaction becomes important, although it is now generally accepted that viscous-inviscid inter- 
action does not suppress the ultimate boundary-layer singularity leading to vortex eruption. However, the inter- 
active stage may certainly be one segment of the vortex eruption and could be used as a means of controlling the 
development of the singularity. This, of course, assumes that marginal states could be found and that some form 
of artificially induced interaction could be used to control those marginal states (a suggestion which seems 
reasonable given recent computations). It seems clear that computational methods may have to accommodate 
the viscous-inviscid interaction at some point in the vortex eruption and such interactions may include artificially 
induced interactions if an attempt is made to control the development of the unsteady flow. It is important, 
therefore, to understand any factors which may complicate the theory and computation of these interactive 
structures - which brings us to the crux of the present study. 

In 1985, Smith Sc Burggraf showed that two-dimensional Tollmien-Schlichting waves can pass into a high 
frequency regime with elevated amplitudes in which the wave is predominantly inviscid. Viscous effects are 
confined to a sublayer which is governed by classical unsteady boundary layer equations and is driven by the slip 
velocity of the inviscid sublayer. As such, the viscous sublayer could admit the Van-Dommelen Si Shen (1980) 
singularity and burst into the outer layers (see also Elliott et al (1983)) . However, with this exception, the viscous 
sublayer is decoupled from the rest of the structure. The inviscid sublayer is found to be governed by nonlinear 
thin-layer Euler equations and it is this layer which feels the effect of the pressure-displacement interaction. For 
hypersonic flows, the predominantly inviscid flow is governed by a modified inviscid Burgers equation 

A r + (A - l)Ax = 0 . (1.1) 

This equation also governs the unsteady initial value problem for high frequency, large amplitude, short-scaled 
waves introduced into the original viscous interactive boundary-layer. Additional implications of eqn. (1.1), 
especially as regards finite-time singularities, are addressed by Brotherton-Ratcliffe Si Smith (1987) and Smith 
(1988). 
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From the Smith & Burggraf (1985) model, eqn. (1.1), several interesting observations can be made. Linear 
waves propagate upstream, finite-time shocks are inevitable at finite wave amplitudes, and at sufficiently high 
amplitudes the characteristics change direction. This suggests that the differencing of at least the streamwise 
convective term uu x and pressure gradient, which produce the (■A - l)Ax term, is not a simple matter, but 
depends on an interplay between the deviation from a "mean flow profile” and the pressure displacement 
interaction. 

Most finite difference methods for calculating high Reynolds number unsteady viscous flows tend to rely on 
upwinding schemes based on the sign of the streamwise velocity (see for example Keller (1978)). However, the 
above observations suggest that traditional differencing techniques may not be appropriate in some unsteady 
interactive boundary-layer computations. The present study indicates that in non-parallel, unsteady, hyperson- 
ic boundary layers — and in a wide variety of zero— displacement boundary layers — the zones of influence are 
determined by a subtle interplay between the convective effects, the pressure-displacement interaction and the 
nonparallel (possibly separated) flow - throughout the entire viscous sublayer. 

2. Governing Equations 

The principal problem to be addressed in this study is the linear unsteady flow superimposed upon an 
originally nonlinear steady hypersonic boundary layer. The equations governing a high Reynolds number un- 
steady interactive hypersonic boundary layer may be found in a number of studies, including: Brown et al 
(1974,1975) and Gajjar & Smith (1985). The governing equations, with the Prandtl transposition, are found to 
be: the conservation of mass equation, 

U x +V r = 0 (2.1) 

and the conservation of streamwise momentum equation, 


U T +UU X +VUy= -P^X, T)+Uyy , ( 2 - 2 ) 

with no-slip boundary conditions at the wall and 

U[X, Y, T) — Y + F{X) + A{X, T) as Y -+ ° ° . (2.3) 

The pressure displacement relation in the hypersonic flow of Brown et al (1974,1975) is simply P=-A. Here 
F(X) is a prescribed steady hump/indentation shape. For the unsteady internal flow of Smith (1976) and Duck 
(1979,1985), as well as the short-scale hump in a Blasius boundary layer (Smith et al (1981)), the P=-A law is 
replaced by A=0. The primary issue of interest here is the numerical integration of the above system of equations 
using a finite difference procedure. This issue, as will be shown in this study, may not be as simple as it might first 
appear. 


3. High Frequency Limit 

Suppose that a nonlinear steady state solution has been calculated for equations (2.1) through (2.3) - say a 
hump induced separation. This may be easily done using a variant of the Davis (1984) alternating direction 
explicit (ADE) method (see also Brotherton-Ratcliffe (1987)). A linear unsteady disturbance is superimposed 
upon this steady base solution: 


[U, V] - [U 0 , K 0 ](X, y) + i [«i, v t ,l(X, Y,7)+ 0{l 2 ) [e <Z 1) . (3. 1) 

Equations (2.1) and (2.2) will be further simplified by letting the time scale become short, or 

— =Q — (Q > l) . (3.2) 

dT dt 0 

The reader will note that this approach follows very closely the work of Smith & Burggraf (1985), Tutty & 
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Cowley (1986), Brotherton-Ratcliffe & Smith (1987) and Smith (1988). The lower deck is found to break up 
into two distinct regions (I and III in Fig. 1) . For purposes of discussion, and comparison with Smith & Burggraf 
(1985), an extra region, II, will be introduced. Region II can be derived as a limit of Region III. The entire 
structure is found to be confined to a neighborhood of a point X 0 within the lower deck, with 

X=X 0 + Q- l x 0 , (3.3) 

and viscous effects are found to be confined to a thin classical Stokes layer, Region I, where Y = Q _1/2 r. The 
Stokes layer equations are solved subject to the no-slip boundary conditions at the wall and matching with 
Region III. It should be noted that Region I is decoupled from Region III. In Region III we take Y to be 0(1) and 
the governing equations are found to be 

“1^0 + = 0 and Wi, 0 + c/ o(y)M lxo + I/'o(y)vi = -Pi X0 , (3.4a, b) 

where U 0 (Y) is the local steady velocity profile. Upon integration of equation (3.4b), application of the tangency 
boundary condition, and the hypersonic pressure displacement relation, the entire problem is reduced to: 

Y 

« lf0 + Uo(Y)u Uq - l/'oM \ u lxQ (xo, * to)dr, = d ljt() _ (3.5) 

0 

subject to U\ Q\ as Y —*■ oo . Equation (3.5) will form the basis of the numerical calculations presented later in 
this study. It may be shown that in the limit as Y becomes large eqn. (3.5) reduces to the linear version of the 
Smith & Burggraf (1985) wave equation (1.1), but with a non-trivial displacement effect from region III: 

d 1/o + (Fo + A 0 -l)d ljro = -(5 . (3.6) 

Further details of this structure may be found in Tutty & Cowley (1986) and Rothmayer (1990). 

4. Zones of Influence for the Unsteady Linear and Nonlinear Viscous Sublayer 

Although the numerical results of this study will be for the high frequency limit of section 3, this approach 
may be generalized to other cases. First consider an integrated form of equations (2.1) and (2.2) for the fully 
nonlinear problem: 

y 

U T + UU X - U Y J UJ,X, r}, l)drj = -PjjC, T) + Uyy , (4.1) 

0 

where P(X,T) is fixed from the appropriate pressure displacement law and the u-matching condition (2.3). 
Linearized problems may be considered, in which case, given an expansion of the form of eqn. (3. 1), equation 
(4.1) becomes 

y 

u T + UpUx + uUox + VqUy - UpyJ UxjX, J], T)drj = -p^X, T) + Uyy , (4.2) 

0 

subject to no-slip at the wall and u(X, oo, 7) = a 4- / as Y -> °o. The high frequency limit of Rothmayer (1990) 
results in a parallel flow approximation and neglecting viscous effects in the main portion of the boundary layer. 
Equivalently, the high frequency equation is simply equation (4.2) without the underlined terms, for linear 
flows, or equation (4.1) without the underlined terms, for nonlinear flows. Equation (4.2) may be evaluated 
on the Y grid Yj — (/-1)AY using, say, a central difference approximations for u Y and Uyy and a trapezoidal 
rule quadrature for the integral in (4.2). As will be shown in section 5, equation (4.2) may be reduced to 
the form 


u r + Au* = Bu + g or [u 1 ], 0 + A[u,] X() = 0 (4.3a,b) 

when a high frequency approximation is used (see eqn. (3.5), also Rothmayer (1990)) and the underlined 
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terms in equation (4.2) are neglected. This equation assumes that a suitable pressure displacement relation 
has been used to relate the pressure to some combination of the velocities on the vertical grid (see Rothmayer 
(1990)). For example, the hypersonic interaction law gives P\ = - U\ iN whereas the prescribed displacement 
law (A=0) may be used, in which case the high frequency approximation gives: 

N-\ 

Px = ~g + ^ and ${x>t) = f + (Y N + A 0 + F 0 - l/2]f x . (4.4) 

2 


The vector g is given by g — [g...g] r while the matrix A is given in Rothmayer (1990) for the high frequency 

limit. The same analysis may be applied to equation (4.1), but now in terms of U = [U iu ...,£/i*] T and the 
coefficient matrix A is a function of the solution. Equation (4.3a) is simply an N-dimensional wave equation 
which may be solved by standard means using the method of characteristics. Note that the viscous effect as 
well as the nonparallel effects (i.e. all underlined terms in (4.2)) only contribute to B and so do not affect 
the sub-characteristic analysis. The primary difference between the present study and the standard method 
of sub-characteristics is that the characteristics cannot be determined via a local analysis, due to the pressure 
displacement interaction P=-A (or the A=0 law). The eigenvalues of the coefficient matrix A can readily be 

calculated and satisfy |A — A^lj — 0, and the eigenvectors are found from Av/ — Xpfi. The actual method used 
in this study for finding the eigenvalues and eigenvectors is a nonlinear Newton—Raphson method (which is 
discussed in Rothmayer (1990)). It is not assumed that the eigenvectors have been normalized. Using classical 

methods (see John (1971)), a new solution vector u = Vu is defined, where V is an NxN matrix whose columns 
are the eigenvectors, i.e. V = [>q ...V/y] . Substitution into equation (4.3a) and multiplication by the inverse 
of V diagonalizes eqn. (4.3a), and gives 

u r + Au* — Bu + g , (4-5) 

where B = V _1 BV, g = V' ] g and A is the diagonal eigenvalue matrix A = diag{X „ Equation (4.5) can 

be easily integrated along its characteristics for u , and then the u can be converted back to u. Direct integration 
of (4.3a) requires a differencing formulation which correctly incorporates the diagonalization transformation 

u - Vu. 


5. A Solution of the Non-Local High-Frequency Characteristic Problem 


As in the preceding section, consider a quadrature of the integral in equation (3.5) on the vertical grid 

Yj = (/-1)AY with j=l,...,N. A trapezoidal rule wall be employed here, although other quadratures may be used, 

in which case: 


“ \ ~ (y AyX s) An/ °/% o + ( U °j 2 An/ °i) 




= 0 


(5.1) 


keeping in mind that = d \ . Equation (5.1) is applied at the points j— 3,...,N. At j-2 the summation from 
L=2,...,j-1 is removed from the equation, while at j=l equation (5.1) is replaced by iij -e ^i Nj(q = 0. The 

above equations are just the single matrix w r ave-equation (4.3b) where [u i] = [wi 1 ... Wi^] an ^ coeffi- 
cient matrix A is given in Rothmayer (1990). The characteristic slopes, Xj, are the eigenvalues of A, and satisfy 
|A - A;l| = 0, where X = dx 0 /dt 0 • The Riemann invariant Tj associated with the A, eigenvalue satisfies (see Cho- 
rin & Marsden (1979)): 
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In general there are N eigenvalues, and it turns out that the first N-l are roughly the slope of the U 0 velocity at a 
particular j gridpoint (for the particular solution being calculated in this study) . The last eigenvalue, does not 
appear to have any simply determined value and must be found via numerical computation. These eigenvalues 
have associated with them the Tj Riemann invariants, each of which is constant along the characteristic with 
slope Xj. Therefore at each (xo,t 0f Y) gridpoint a system of N equations in the N u } . ’ s needs to be solved, namely 


z 



(5.3) 


where (r ,) 0 ’s are the Riemann invariants evaluated along the initial data plane. The equations (5.3) may be 
inverted to give: 


m . (5.4) 


This is effectively the diagonalization process of equation (4.5). Equation (5.4) gives the dependence of each 
uij on the Riemann invariant T ; associated with the eigenvalue. Therefore the inverse of the Jacobian matrix 

of the Riemann invariants gives the domain of dependence of the streamwise velocity at each jth gridpoint in 
terms of the j characteristics with slope Xj t providing that the Riemann invariant multiplying a given row element 

of the inverse is non-zero. In the general case, there is no reason to expect the domains of dependence to be 
simple. 


For purposes of computation, an idealized flow will be considered. The velocity field is assumed to take the 
form: 


U 0 {Y) = V+A(l-e- i ) t 


(5.5) 


where A is taken as an independent parameter, in lieu of A 0 +F. The wall shear stress is given by 
*w = Uq[0) = 1 + A, indicating that the flow is attached for A>-1 and separated for A<-1. The case A=0 gives an 
undisturbed sublayer. Velocity profiles for various values of A are shown in Figs. 2 through 4, along with the 
eigenvalue and eigenvector calculations. In addition, a new function, 0,, will be defined, which is the product of 
a row element in the inverse eigenvector matrix and the corresponding element of the Riemann invariant vector: 


0, - 


ar 

du { 


Mo 


(5.6) 


The results of these calculations are shown in Figs. 5 through 8, assuming constant perturbation velocities on 
the initial data plane. The results for A=0 are in agreement with the linear results of the Smith & Burggraf (1985) 
study. The solution at any vertical point in the grid depends only on the j=N eigenvalue which has slope -1. The 
results of the above calculations for an accelerated flow are shown in Figs. 6 and 7. Consider the schematic 
interpretation of Figs. 5 through 7, shown in Fig. 8. This figure is a qualitative interpretation only and is not 
meant to convey accurate quantitative data. At low amplitudes (i.e. A near 0) all of the characteristics point 
downstream with slope -1. As the flow is accelerated (i.e. A increasing) a region of dependence begins to 
emerge for small values of j, and the slope of the downstream-directed characteristic begins to decrease. A 


233 


typical case, say A=5, now has the j=N eigenvalue competing with the eigenvalues which are approximately the 
velocities near the bottom of the boundary layer. Note that this only occurs in the outer portion of the boundary 
layer (i.e. for j large) . The eigenvalues associated with small j’s have characteristics which point upstream and so 
the tendency in the outer portion of the boundary layer is for the upstream directed characteristics to begin 
competing with, and eventually overtaking, the downstream directed characteristic. The overall picture is that 
the characteristics in the outer portion of the sublayer appear to change direction for increasing A, whereas the 
characteristics in the lower portion of the sublayer do not. This means that the characteristics are not pointing in 
the same direction throughout the sublayer when the change in the direction of the characteristics does occur, 
but vary with vertical position. 

7. Conclusion and Implications for Finite Difference Computations 

In this study it has been shown that the sub-characteristics in an unsteady interactive viscous flow are not 
simply determined by the convective velocity, but rather are fixed by an interplay between the convective terms, 
the pressure displacement interaction, and the nonparallel base flow. In addition, the sub— characteristics may 
vary throughout the entire viscous sublayer and may possess a complex structure. These results seem to be in 
accord with the work of Smith & Burggraf (1985) on nonlinear hypersonic waves, which suggests that the 
characteristics will change direction at sufficiently large disturbance amplitude. Is it worthwhile to perform the 
domain of dependence calculation before proceeding to the finite difference solution? This of course depends 
on the results for a particular case. A simple differencing scheme may correctly capture the physics of the 
problem in question. How'ever, the present study indicates that it is possible for an unsteady flow to possess 
complicated domains of dependence, and hence to require complex differencing schemes. It is anticipated that 
problems related to incorrectly modeling the domains of dependence in a finite-difference method will manifest 
themselves either as a CFL (Courant-Friedrichs-Lew'y) condition or as spurious oscillations in the solution. 
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Figure 1. Asymptotic structure for a high frequency 
unsteady linear disturbance in a nonparallel viscous 
sublayer. 
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Figure 2. Eigenvalues, eigenvectors and the velocity 
profile for A=0 and N=51. ej is the eigenvector 
associated with X;. 


O : eigenvalue A ; 
— : velocity Uq(Y,) 
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Figure 3. Eigenvalues, eigenvectors and the velocity 
profile for a model separated flow (A=-5, N=51). 


Figure 4. Eigenvalues, eigenvectors and the velocity 
profile for a model accelerated flow (A=5, N=51). 
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Figure 5. Domains of dependence for vertical 
positions in the grid and A=0. 


Figure 6. Domains of dependence for various vertical 
positions in the grid and A=5. 
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Figure 7. Domains of dependence for various vertical 
positions in the grid and A=10. 


Figure 8. Schematic diagram of the domains of 
dependence given by Figs. 5 through 7. Results are 
qualitative only. 









